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A nonlinear Landau-Zener model was proposed recently 
to describe, among a number of applications, the nonadia- 
batic transition of a Bose-Einstein condensate between Bloch 
bands. Numerical analysis revealed a striking phenomenon 
that tunneling occurs even in the adiabatic limit as the non- 
linear parameter C is above a critical value equal to the gap 
V of avoided crossing of the two levels. In this paper, we 
present analytical results that give quantitative account of 
the breakdown of adiabaticity by mapping this quantum non- 
linear model into a classical Josephson Hamiltonian. In the 
critical region, we find a power-law scaling of the nonadi- 
abatic transition probability as a function of C/V — 1 and 
a, the crossing rate of the energy levels. In the subcritical 
regime, the transition probability still follows an exponential 
law but with the exponent changed by the nonlinear effect. 
For C/V » 1, we find a near unit probability for the transi- 
tion between the adiabatic levels for all values of the crossing 
rate. 

PACS: 32.80.Pj, 03.75.Fi, 73.40.Gk, 03.65.-w 



I. INTRODUCTION 

It is a common practice in the study of quantum sys- 
tems to consider only a finite number of energy levels 
which are strongly coupled. The special case of trunca- 
tion to two relevant levels is of enormous practical in- 
terest, and a vast amount of literature has been devoted 
to the dynamical properties of the two-level systems [0. 
The Landau-Zener tunneling between energy levels is a 
basic physical process Q, and has wide applications in 
various systems, such as current driven Josephson junc- 
tions atoms in accelerating optical lattices and 
field-driven superlattices 

A nonlinear two-level system may arise in a mean field 
treatment of a many-body system where the particles 
predominantly occupy two energy levels, where the level 
energies depend on the occupation of the levels due to 
the interactions between the particles. Such a model 
arises in the study of the motion of a small polaron 
for a Bose-Einstein condensate in a double- well po- 
tential 10-^ or in an optical lattice |10-12|, or for a 
small capacitance Joseph-junction where the charging en- 
ergy may be important. In contrast to the linear case, 
the dynamical property of a nonlinear two-level model 
is far from fully understood, and many novel features 
are revealed recently p3|,H , including the discovery of a 



nonzero Landau-Zener tunneling probability even in the 
adiabatic limit when the nonlinear parameter exceeds a 
critical value. 

In this paper, we present a comprehensive theoretical 
analysis of the nonlinear Landau-Zener tunneling, ob- 
taining analytical results for the tunneling behavior in 
various regimes. This is made possible by the observation 
that the population difference and the relative phase of 
the two levels form a pair of canonically conjugate vari- 
ables of a classical Hamiltonian, the Josephson Hamil- 
tonian 1^,^. The fixed points of this classical system 
correspond to the eigenstates of the nonlinear two-level 
Hamiltonian. Adiabatic evolution of the fixed points as 
a function of level bias correspond to adiabatic evolution 
of the two eigenstates. The breakdown of adiabaticity 
occurs when a fixed point turns into a homoclinic orbit 
which is possible if and only if the nonlinear parameter 
exceeds a critical value. Then the transition probabil- 
ity between the two energy levels can be nonzero even 
in the adiabatic limit, and is found to rise as a power 
law of the nonlinear parameter minus its critical value. 
Right at the critical value of the nonlinear parameter, the 
transition probability goes to zero as a power law of the 
crossing rate of the two levels. Below the critical point, 
the transition probability follows an exponential law as 
in the linear case but with the exponent changed due to 
the nonlinearity. Far above the critical point, we find a 
near unit probability of transition between the adiabatic 
levels for all values of the crossing rates. 

Our paper is organized as follows. In Sec. II we in- 
troduce the nonlinear two-level model and discuss the 
behavior of its eigenstates and eigenenergies. In Sec. HI, 
we cast the nonlinear two-level system into the Josephson 
Hamiltonian and study the evolution of the fixed points 
as functions of the system parameters. In Sec. IV we 
derive a formula for the tunneling probability in the adi- 
abatic limit and analyze its power law behavior above 
and near the critical value of the nonlinear parameter. 
In Sec. V, we calculate the nonadiabatic tunneling prob- 
ability in the critical and subcritical regimes. In Sec. VI, 
we derive the nonadiabatic tunneling probability for the 
ultrastrong nonlinear coupling using the stationary phase 
approximation. In sec. VII, we draw our conclusions and 
discuss how our findings may be observed experimentally. 
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II. THE NONLINEAR TWO-LEVEL MODEL 

Our model system consists of two levels as in the stan- 
dard Zener model but with an additional energy differ- 
ence depending on the population in the levels, ft is 
described by the following Hamiltonian matrix 



i\b\' 

V 
2 



V 



(1) 



where a and b are the probability amplitudes. The Hamil- 
tonian is characterized by three constants: the coupling 
between the two levels V, the level bias 7 as in the linear 
Zener model, and the nonlinear parameter C describing 
the level energy dependence on the populations. The 
amplitudes satisfy the Schrodinger equation 



(2) 



\b\' 



according to which, the total probability |a| 
conserved and is set to be 1. 

We wish to study how the system evolves when the 
level bias 7 changes slowly from —00 to -l-oo. It will 
be useful to find the adiabatic levels by diagonalizing 
the Hamiltonian (]^), i.e., solving the following eigen- 
equations with an eigenenergy e, 



,1 



\af))a + ^b 



ea. 



(3) 



(4) 



For a solution with nonzero amplitudes, we impose the 
determincntal condition 



det 



V 

2 



V 



\an-e 



(5) 



which leads to a quartic algebraic equation for the 
eigenenergy. 



1 T 



4 



16 



0. 



(6) 



The eigenenergies e correspond to the real roots of the 
above quartic equation. It is readily found that this quar- 
tic equation has two real roots when C < V, while there 
are four real roots when C > V . Two typical param- 
eters are chosen to demonstrate the two cases in Fig.l. 
At C /V — 0.5(Fig.la), there are two energy levels corre- 
sponding to the two real roots; At C /V = 2(Fig.lb), a 
loop appears at the tip of the lower level in the regime 



—Ic < 7 < 7c, where Eq. (g) has four real roots (this can 
be easily verified from the above equation at 7 = 0). The 
relation between the singular point 7c and the parame- 
ters will be given by Eq. (|l9|) in the following section. 

The eigenstates can be determined at the same time 
from the above equations. 



7 



2C + 4e' 



6 = 



V^/2CTTt 
4eVC-' + 2e + 



(7) 



(8) 



These eigenstates are not orthogonal to each other for 
finite 7, but become so in the limits of 7 — + ±cxd, where 
e ±|7|/2. For instance, for the lower level, we have 
(a, 6) (1,0) at 7 — > —00 and (a, 5) — > (0,1) at 7 ^ 
+00. Due to the nonlinearity of the model, there are two 
more eigenstates for C > and in the region —7c < 
7 < 7c. These correspond to the loop structure in the 
adiabatic energy levels shown in Fig. lb. 



III. THE CLASSICAL ANALOG 

In this section, we cast the nonlinear two-level sys- 
tem into an effective classical system, i.e. the Joseph- 
son Hamiltonian. The dynamical mechanism for the 
breakdown of the adiabatic evolution will be revealed 
by investigating this classical Hamiltonian system. The 
Schrodinger equation (0) has the following explicit form 



(9) 



'1 
'2 



\a\'))b- 



V 



(10) 



'^'which describes the evolution of the amplitudes for gen- 
eral occupations and phase of the wavefunctions. The 
physically interesting information of the state is com- 
pletely represented by the population difference 



s{t) = \b\'^\a\^ 
and the relative phase 



(11) 



(12) 



From Eqs. (|9|) and (10), we can obtain the following 
equations of motion 



dt 



s{t) = -Vy/l-s^{t) sin6'(i). 



Cs{t) 



Vs{t) 



COS 9{t). 



(13) 



(14) 
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These equations can be rewritten in the canonical form 
s 

nian given by 



- ^mf- , = , with an effective classical Hamilto- 



C 



He{s,d,7) = + 7s(t) - V^l - S^t) COS0{t), 

(15) 

which has the form of a Josephson Hamiltonian. 

Since we are interested in the situation of a quasistatic 
level bias 7, it will be useful to find the fixed points of 
the classical dynamics when 7 is fixed in time. The fixed 
points correspond to the eigenstates of the nonlinear two- 
level system, and are obtained by equating the right hand 
sides of eqs.(O) and (nj) to zero, yielding 



= 0, t:, 



7 + Cs* + 



Vs* 



■.COS0* = 0, 



(16) 



(17) 



VI -s*' 

These can be consolidated into a single quartic equation 

72 



s*' + 2ls*' - (1 - ^ - - 2ls* ~^ = 0. 



(18) 



The number of the fixed points depends on the nonlin- 
ear parameters C. For weak nonlinearity, C/V < 1, only 
two fixed points exist, denoted as Pi and P2 in the pan- 
els of Fig. 2. These are elliptic fixed points corresponding 
to the maximum and minimum of the classical Hamil- 
tonian, and each arc surrounded by closed (elliptic) or- 
bits. The fixed points are located on the lines of 6** = tt 
and 0, meaning that the two corresponding eigenstates 
of the two level system have tt and relative phases. As 
the level bias changes from 7 = —00 to +00, Pi moves 
smoothly along the line 6* — n from the bottom (s = —1) 
to the top (s = -l-l), corresponding to the lower branch 
of the energy levels in Fig. 1(a). On the other hand, P2 
moves from the top to the bottom corresponding to the 
upper branch. 

For stronger nonlinearity, C/V > 1, two more fixed 
points can appear when the level bias lies in a window 
~7c < 7 < 7c- The boundary of the window can be ob- 
tained from the condition of real roots for quartic equa- 
tions, yielding 



7, = (C2/3_V^2/3)3/2^ 



(19) 



As shown in Fig.3(c-e), both of the new fixed points lie 
on the line 9* = tt, one being elliptic (P4) as the origi- 
nal two and one being hyperbolic (P3) corresponding to 
a saddle point of the classical Hamiltonian. One of the 
original fixed point, P2, still moves smoothly with 7, cor- 
responding to the upper adiabatic level in Fig. 1(b). The 
other, Pi, moves smoothly up to 7 = 7c, where it collides 



with P3, corresponding to the branch OXT of the lower 
level in Fig. 1(b). The new elliptic point P4, created at 
7 = —7c together with P3, moves up to the top, corre- 
sponding to the branch WXM of the lower level. The 
hyperbolic point P3, moves down away from its partner 
after creation and is annihilated again with Pi at 7 = 7c, 
corresponding to the top branch WT of the lower level. 
The eventual fate of the fixed point Pi will determine the 
adiabatic tunneling probability as shown below. 



IV. ADIABATIC TUNNELING DUE TO 
NONLINEARITY 

For quasistatic change of the level bias 7, a closed orbit 
in the classical dynamics remains closed such that the 
action 



(20) 



is invariant in time according to the classical adiabatic 
theorem [|l^ , which is valid as long as the relative change 
of the system parameter in a period of the orbit is small. 



Tdj/dt < 7. 



(21) 



The action equals the phase space area enclosed by the 
closed orbit, and is zero for a fixed point which has no 
area. Since the closed orbits surrounding an elliptic fixed 
point all have finite periods T, they should evolve adi- 
abatically with the area of each fixed in time. We thus 
expect an elliptic fixed point to remain as a fixed point 
during the quasistatic change of the system parameter. 
For the case oi C/V < 1, the two fixed points (both ellip- 
tic) evolve adiabatically throughout the entire course of 
change in the level bias, implying the absence of transi- 
tion between the eigenstates in the adiabatic limit. This 
is still true for the fixed point P2 in the case of C/V > 1, 
meaning a state starting from the upper level will remain 
in the upper level. 

The adiabatic condition is broken, however, when Pi 
collides with the hyperbolic fixed point P3 to form a ho- 
moclinic orbit where the period T diverges. Nevertheless, 
the classical 'particle' will remain on this orbit, because 
the orbit is surrounded from both outside and inside by 
closed orbits of finite periods, which should form adia- 
batic barriers to prevent the particle from escaping. Af- 
ter this crisis, the homoclinic orbit turns into an ordinary 
closed orbit of finite period, and will evolve adiabatically 
for 7 > 7c according to the rule of constant action which 
is now nonzero. In this region of the level bias, the pop- 
ulation difference oscillates along the orbit, which is just 
the feature described in Ref. that "there is a violent 
shaking above and about the lower branch of the adia- 
batic level after the terminal point" . 



3 



The adiabatic tunneling probability can then be cal- 
culated from the action Ic of the homoclinic orbit. This 
orbit eventually evolves into a straight line at s = s/ (see 
Fig. 4), where the condition of constant action demands 
that 



s/ = l-/(s,). 



(22) 



On the other hand, from the definition of population dif- 
ference, we find the final probabilities to be given by 



\bf? 



(23) 



The adiabatic tunneling probability is then expressed as 

r,d = |a/|2=/,/2. (24) 

The task we face then is to calculate the action of the 
homoclinic orbit. Let Sc be the population difference at 
the degenerate point Pc where Pi collides with P3 when 
7 = 7c. It can be obtained by substituting the expression 
(|T|) into Eq. (|T|) as 



(25) 



Considering the fact that the degenerate point Pc lies on 
9 = TT, we find the total energy (the value of the classical 
Hamiltonian) as 



C 



(Sc)^ + 7cSc + Vy/1 - (Sc)2 



(26) 



The homoclinic orbit should have this energy, and its 
trajectory 



s = s{9;Ec), 



(27) 



may be obtained by equating the classical Hamiltonian 
( llSl ) to Ec- The corresponding action can then be evalu- 
ated to give the tunneling probability as 



(28) 



This has been evaluated numerically, and the results com- 
pare extremely well with those obtained by directly inte- 
grating the time-dependent nonlinear Schrodinger equa- 
tion using the Runge-Kutta algorithm (See Fig. 5). 

The adiabatic tunneling probability can be evaluated 
analytically in the critical region oi 5 = C /V —1^0. 
The singular point of the level bias is found to leading 
order as 



7c 



(29) 



The homoclinic orbit is confined near the critical point, 
with its top at 



St 




(30) 



We may expand the classical Hamiltonian to leading or- 
ders of s — Sr and 9 — it to find 



1 



(31) 



From the area of this orbit the adiabatic tunneling prob- 
ability for this limiting case is found to be given by the 
power law 



f^' 4 3 

Tad^ / {9-n)ds = —5-^ 
Jo 37r 



(32) 



We may also obtain the adiabatic tunneling probability 
analytically in another limit, C/V — > 00. To leading 
order in the small quantity a = V/C, we find 



7e C(l - '-a'/% 



(33) 



(34) 



In this case, the homoclinic orbit is confined near s = — 1, 
and the deviation is the same order as Sc — (~1)- If we 
write s — (—1) = rf'a^l^ , then 



3 +77^ - 2 \/2 r;cos(6i) ~ 0. 



(35) 



By solving this equation and considering s = ■r]^(9)a'^l'^ — 
1, the adiabatic tunneling probability can be obtained 



lad =i 1 - ■ 



(36) 



V. NONADIABATIC TUNNELING AT AND 
BELOW THE CRITICAL POINT 

We have thus found that it is possible to tunnel be- 
tween the adiabatic levels even in the adiabatic limit 
if the nonlinear parameter exceeds a critical value, i.e., 
C jV > 1. We now consider the case of nonzero sweeping 
rates, and study how the nonliearity affects the nonadia- 
batic tunneling probability. In the linear case C = 0, the 
Landau-Zener formula prescribes an exponential tunnel- 
ing probability between the adiabatic levels. 



r 



exp(. 



"2^ 



(37) 



For the nonadiabatic case, we will focus our attention 
in this section to the critical point and the subcritical 
region, and consider the near adiabatic case (i.e., a 7^ 
and a « 1). 
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For this purpose, we need to investigate the evolution 
of the fixed point Pi as well as the nearby periodic orbits. 
Then, in addition to the action variable / mentioned in 
the above section, we should introduce its canonical con- 
jugate - the angle variable (f>. They satisfy the following 
differential equations p5[ |. 



dR 



7, 



1 rr , 9R. 

4' = ^(-^,7) + -gjl 



(38) 



(39) 



The function R{I,4>) is defined as J s{9,E,j)d9, 
where E is the energy of the periodic orbit. Here, uj 
is the frequency of orbits in the immediate neighborhood 
of the fixed point Pi. Its expression is particular sim- 
ple at / = 0, corresponding to the fixed points Pi. It 
can be calculated by linearizing the equations of motion 
eqs.( [r^JT^ ) near the fixed point as 

cu* ^ V^l/{l-is*f)~^^l-{s*f. (40) 

As in the adiabatic case considered in the last section, 
the transition probability is still given by the increment 
of the action, i.e., 



r = -Ai. 

2 



(41) 



The reason is as follows. The initial state at 7 = —00 is 
the fixed point Pi with zero action, and the final state is 
an orbit of finite action which becomes a horizontal line 
s = at 7 = -|-oo. Then an integration of Eq. (^8|) yields 



AI 



+°° dRdrf_, 
86 dt 



(42) 



To evaluate this integral, we need to express as a 
function of (j) itself. Near the adiabatic limit, we may 
omit the second term in Eq. (|39|) and set uj{1 , 7) = i^*(7) 
in the first term, yielding 



(43) 



On the other hand, by substituting 9* = tt into eq.(|T^) 
and differentiating with respect to time, one gets 



dt 
ds* 



(44) 



Combining these equations, one may relate s* to (j) a-nd 
thus express as a function of </> itself. 

The principal contribution to the integral comes from 
the neighborhood of the singularities of the integrand, 



which comes only from the zeros of the frequency cj) = 
uj*{'y). These zero points are easily found from Eq.(|40|) 
as 



(45) 



As will be shown below, the integral is exponentially 
small if there is no real singularities, and becomes a power 
law in the sweeping rate if there is a singularity on the 
real axis. 

We first consider the case of critical nonlinearly, 
C/V — 1, for which the singular point occurs at s* = 0. 

Near this point, we find from Eqs.(^4|) that uj* z 

and (j) ~ g:(|)'^^^-^(s*)^- Then, we have an approximate 
relation uj* ^ ai(j)i near the singularity. Substituting 
these expressions back to the eq.(^2|), and utilizing the 
fact that ^ is independent of a, we find a power-law 
behavior for the tunneling probability 



IVs* 



T ■ 



3 



(46) 



This scaling law has been verified by our numerical cal- 
culations (Fig. 6). 

Now we extend our discussion to nonadiabatic tunnel- 
ing for subcritical nonlincarity, where the zeros of the 
frequency oj* are complex. Because i? is a periodic func- 
tion of the angle variable, it can be expanded as a Fourier 
series, 



dR 



+ CXD 



^ ile'^'f'Ri. 



(47) 



After substituting this series into eq. (^ , we may deform 
the contour of integration into the complex plane. For 
those positive / terms in the series, we may raise the 
contour into the upper half plane until it is "caught up" 
by a singularity of the integrand, namely a zero oiui* . Let 
00 be the singularity nearest to the real axis, i.e. the one 
with the smallest positive imaginary part. The principal 
contribution to the integral comes from the neighborhood 
of this point, yielding [hsl. 



AI - exp(-;jm((?;>o)), 



(48) 



where 'Im((/)o)' indicates the imaginary part of the com- 
plex number and Is is the lowest order of the Fourier 
series whose Fourier component does not vanish. The 
negative I terms give the same exponential and so affects 
only the prefactor which we do not consider here. In the 
neighborhood of the elliptic fixed point Pi, the system 
can be approximated by a harmonic oscillation with the 
function R ^ sin(2(/)), so we have Is = 2. 

To find the complex zero (/iq , we combine the equations 
of (H) to express it as an integral of s* 



a 



5 



{s*ff'^~C/V\ds 



(49) 



We take the upper bound to be the zero Sq of the fre- 
quency, while the lower bound only reflects a certain 
choice of the origin of (f) which does not affect the imagi- 
nary part. 

The tunneling probability is thus found to be propor- 
tional to the exponential 



problem becomes adequate. We begin with the variable 
transformation 

a = a'exp(-z j\l + ^{^^^ - \a\^))dt) , (52) 



V ~ exp(-g-^), 
where the factor in the exponent is given by 



(50) 



q=^ {1 + xy/\l/{l + x^f/' - C/Vf/'dx 
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(51) 

For the linear case C = 0, the factor q is exactly unit, 
which means that the expression of the tunneling proba- 
bility is identical to the standard Landau-Zener formula. 
For the nonlinear case, C/V > 0, this factor becomes 
smaller than one, indicating the enhancement effect on 
the nonadiabatic tunneling. As C/V goes up to 1, the 
critical point, this factor vanishes, meaning the break- 
down of the exponential law. 

In Fig. 7, we plot this factor and compare it with 
the results of numerical integration of the nonlinear 
Schrodinger equation. The agreement is very well for 
the nearly linear case, C/V << 1, and there is some 
deviation when C/V is close to unity. One reason for 
the latter is that the prefactor may be important in such 
a case. Also, numerical accuracy may be blamed. To 
obtain the correct slope in the plot of InT vs. a, one 
should calculate the tunneling probability at infinitesimal 
a. However, in practical calculations, this parameter is 
chosen as a finite value instead of infinitesimal. In [ p^ , 
the slope is read in the regime a G [0.02,0.01]. In our 
calculations, this range is extended to [0.005, 0.0025]. As 
is well known, in the regime of small a, the numerical 
fluctuation is serious. This may result in a relative large 
slope numerically. 



VI. TUNNELING AT STRONG NONLINEARITY 

In an earlier section, we have found that for strong 
nonlinearity C/V » 1, there is a near unity tunneling 
probability to the upper adiabatic level even in the adi- 
abatic limit. This probability can only gets larger when 
the sweeping rate is finite. We thus expect that the am- 
plitude b in the original Schrodinger equation remains 
small all the times, and a perturbative treatment of the 



/■* ' 

b = b' exp {i / (■; 

"'0 ' 



C 



(53) 



Then the diagonal terms in Eqs. (y,^^ are transformed 
away, and our above discussions imply that we can take 
a' w 1, and 

rt rt 

dtexp{-i {-f + C{\b\^-\a\^))dt). (54) 



b'^^L 

2i 



We need to evaluate the above integral self- 
consistently. Because of the large C, the nonlinear term 
in the exponent generally gives a rapid phase change 
which makes the integral small. The dominant contri- 
bution comes from the stationary point to of the phase 
where (—7 + C(l - 2|6|^))fo = 0. Expanding about this 
point, we may write 

-7 + C(l-2|6|2) =-«(t_to), 



with 

We thus have 
\bf 



a = a + 2C[^Jb\%^. 



d_ 

'di' 



^ 2 ' 



dtexp{--a{t - tof) 

1 ^ 



(55) 
(56) 

(57) 



We can differentiate this expression and evaluate its re- 
sult at time to, obtaining a few standard Fresnel inte- 
grals with the result [4|6|^]to = Q^Y \f\- Combining 
this with the relation ( pq ) , we come to a closed equation 
for a. 



(58) 



The nonadiabatic transition probability F is given by 



r = i-|5|^ 



^ 2 ' 



1 



dtexp(--a(t - tof) 



2a 



(59) 



Then the above result yields a closed equation for the F, 



1 



1 -F 



1 V2C r- 
P + -V^ 



(60) 



where P = In the adiabatic limit, i.e. 4 = 0, we 



2a 

find that F = 1 — 1.7(^)3. This gives the same exponent 
as the exact asymptotic result (^) obtained in an earlier 
section, but the coefficient differs by about 10%. In the 
sudden limit, — s- cx), we have F — 1 — P, which is exact. 
In figure 8, we compare the above analytical results with 
that from directly solving the Schrodinger equation and 
show a good agreement. 



6 



VII. CONCLUSIONS AND DISCUSSIONS 

In this paper, we present a theoretical analysis on the 
Landau-Zener tunneling in a nonlinear two-level system. 
Our results can be summarized as: 1) The adiabatic 
tunneling probability depends only on the ratio between 
the nonlinear parameter C and the coupling constant V . 
With increasing the ratio C /V , a transition from zero 
adiabatic tunneling probability to nonzero adiabatic tun- 
neling probability occurs at a critical point C /V = 1. 
The critical behavior near the transition point follows 
a power-law scaling. 2) In the subcritical regime, the 
nonadiabatic tunneling probability follows a exponential 
law, which modifies the Landau-Zener formula by adding 
a factor in the exponent. This exponential law breaks 
down at the critical point, instead, a power-law shows 
up. 3)Far above the critical point, we find a near unit 
probability of transition between the adiabatic levels for 
all values of the crossing rates. The analytical expression 
is obtained using the stationary phase approximation. 

The experimental observation of these findings is a 
topic of great interest. It has been shown that the non- 
linear two-level model can be used to understand the 
tunneling of a Bose-Einstein condensate between Bloch 
bands in an optical lattice Another physical model 
for direct application of our theory is a Bose-Einstein 
condensate in a double- well potential [0J|,|l^. The am- 
plitudes of general occupations Ni^2{t) and phases 9i^2 
obey the nonlinear two-mode Schrodinger equations sim- 
ilar to (^) (e.g. see After introducing the new vari- 
ables s{t) = {N2it)~Ni{t))/NT and 6* = 6*2 -6*1, one also 
obtain a Hamiltonian having the same form as Eq.(|l5|) 
except for the parameters replaced by = 2K/h, 7 — 
-[{E° - E^) - [Ui - U2)NT/2]/h, C^{Ui+ U2)NT/2h. 
Here the constants have following meanings: Total num- 
ber of atoms Ni + N2 — Nt, is conserved, 2 are zero- 
point energy in each well, are proportional to 
the atomic self-interaction energy, and K describes the 
amplitude of the tunneling between the condensates. 

In practice, such a double well may be achieved by a 
laser sheet dividing a trap, and the energy levels may be 
moved by shifting the laser sheet. Then, our theoretical 
results can be directly applied to this system without 
intrinsic difficulty. One can design an experiment for 
observing the transition to nonzero adiabatic tunneling 
and verify those scaling law found in our paper. We hope 
our discussions will stimulate the experimental works in 
the direction. 
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Figure Captions 

Fig.l Adiabatic energy levels for two typical nonlinear 
cases (solid line). The dashed lines is the ones for lin- 
ear case(C = 0). (i = 1, • • • , 4) is the fixed points of 
He system, corresponding to the quantum eigenstates as 
shown in (b): OXT Pi, MXW P4,, WT P3, 
only P3 is a unstable saddle point, others are stable el- 
liptic point. 

Fig. 2 The phase-space portrait of the Hamilton sys- 
tem for various parameter 7 at C/V = 0.5. The arrows 
refer to the direction of the shift of the fixed points Pi as 
7 increases. The closed curves are the periodic trajecto- 
ries. In this case, no collision between fixed points occurs, 
which indicates a zero adiabatic tunneling probability. 

Fig. 3 The phase-space portrait of the Hamilton Hf, sys- 
tem for various parameter 7 at C/V = 2. The arrows 
refer to the direction of the shift of the fixed points as 
7 increases. In this case, a collision occurs at a singular 
point 7c, where P2 and P3 collide and form a homoclinic 
orbit with nonzero action. This jump on the action leads 
to a nonzero adiabatic tunneling probability. 

Fig. 4 The homoclinic trajectories associated with the 
singular points Pc at critical point 7 = 7c (left column), 
and its final shape at 7 = +00 after adiabatic evolu- 
tion (right column), respectively. For the action keeps 
as an invariant during the evolution process, the areas 
shadowed are equal for the same C/V. Here we show 
two typical situations corresponding to smaller and larger 
nonlinearity, respectively. 

Fig. 5 The adiabatic tunneling probability vs. param- 
eter C/V. Solid line is our theoretical results, whereas 
the symbols- dot, star and plus denote the results from 
solving nonlinear Schrodinger equation numerically. In- 
set is a magnification of the local plot. 

Fig. 6 The dependence of the tunneling probability on 
the parameter a for the critical nonlinearity. 

Fig. 7 The dependence of the factor q on the ratio C/V. 

Fig. 8 The comparison between our analytical theory 
and the numerical simulations in the regime of strong 
nonlinear coupling. 
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